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1.  Introduction 

The  problem  of  determining  rigid  body  motion  and  surface  structure  from  image  data 
has  been  the  topic  of  many  research  papers  in  the  area  of  machine  vision  [1-22].  Many 
approaches  based  on,  tracking  feature  points  [5,11,19,20]  or  contours  [9],  using  motion 
flow  field  [1,3,4,10,12,16,17,21,22],  texture  [2],  or  image  intensity  gradients  [14,15]  have 
been  proposed  in  the  literature. 

In  the  feature  point  matching  schemes,  information  about  a  finite  number  of  well- 
seperated  points  is  used  to  recover  the  motion  (general  8-point  2-frame  algorithms  of 
Longuet-Higgins  [11],  Tsai  and  Huang  [20],  Buxton  et  al.  [5],  and  the  algorithm  of  Tsai, 
Huang  and  Zhu  [19]  for  planar  surfaces).  These  methods  require  identifying  and  matching 
feature  points  in  a  sequence  of  images.  The  minimum  number  of  points  required  depends 
on  the  number  of  image  frames.  With  2  frames,  in  most  cases,  a  minimum  of  5  points 
results  in  a  unique  solution  from  a  set  of  nonlinear  equations.  However,  using  8  points, 
as  in  algorithms  cited  above,  one  only  solves  linear  equations.  Here,  it  is  assumed  that 
the  more  difficult  problem  of  establishing  point  correspondence  has  already  been  solved. 
In  general,  this  involves  determining  corners  along  contours  using  iterative  searches.  For 
images  of  smooth  objects,  it  is  difficult  to  find  good  features  or  corners. 

For  the  general  case  of  smooth  surfaces,  Longuet-Higgins  and  Prazdny  [11]  suggested 
a  method  that  uses  the  optical  flow  and  its  first  and  second  derivatives  at  a  single  point. 
Later,  Waxman  and  Ullman  [21]  developed  this  into  an  algorithm  for  recovering  the  struc¬ 
ture  and  motion  parameters  from  a  set  of  nonlinear  equations.  Subbarao  and  Waxman 
[17]  recently  found  a  closed  form  solution  to  the  original  formulation  in  [21]  for  planar 
surfaces.  These  methods  while  mathematically  elegant  are  very  sensitive  to  errors  in  the 
optical  flow  data  since  second  order  derivatives  of  noisy  data  are  used. 

At  the  expense  of  more  computation,  more  robust  algorithms  have  been  suggested 
using  the  optical  flow  at  every  image  point  [1,3,4]. 

Longuet-Higgins  [12]  has  presented  a  closed  form  solution  for  planar  surfaces,  very 
similar  to  ours,  using  the  coefficients  of  the  second  order  optical  flow  equations.  However, 
it  is  assumed  that  both  components  of  the  flow  field  have  already  been  computed  for  a 
minimum  of  5  image  points. 

By  representing  a  planar  surface  in  the  form  of  a  closed  contour,  Kanatani  [9]  has 
shown  that  the  surface  and  motion  parameters  can  be  computed  by  measuring  “diam¬ 
eters”  of  the  contour  using  line  and  surface  integrals.  Here,  no  point  correspondence  is 
required. 

Assuming  that  the  planar  surface  has  a  uniform  texture  density,  Aloimonos  and  Chou 
[2]  have  presented  a  procedure  for  computing  the  motion  and  surface  orientation  from 
texture. 

In  much  of  the  research  work  in  recovering  surface  structure  and  motion  from  the 
optical  flow  field,  it  is  assumed  that  a  reasonable  estimate  of  the  full  optical  flow  field 
is  available.  In  general,  the  computation  of  the  local  flow  field  exploits  a  constraint 
equation  between  the  local  intensity  changes  and  the  two  components  of  the  optical  flow. 
However,  this  only  gives  the  component  of  the  flow  in  the  direction  of  the  intensity 
gradient.  To  compute  the  full  flow  field,  one  needs  additional  constraints  such  as  the 
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heuristic  assumption  that  the  flow  field  is  locally  smooth  [7,8j.  This,  in  many  cases,  leads 
to  optical  flow  fields  that  are  not  consistant  with  the  true  motion  field. 

In  an  earlier  paper,  we  presented  an  iterative  scheme  for  recovering  the  motion  of  an 
observer  relative  to  a  planar  surface  directly  from  the  image  brightness  derivatives,  and 
without  the  need  to  compute  the  local  flow  field  [14,15].  Further,  using  a  compact  vector 
notation,  we  showed  that,  at  most,  two  interpretations  are  possible  for  planar  surfaces 
and  derived  the  relationship  between  them.  Here,  we  present  a  closed  form  solution  to  the 
same  problem.  We  first  solve  a  linear  matrix  equation  for  the  elements  of  a  3x3  matrix 
using  intensity  derivatives  at  a  minimum  of  8  non-colinear  points.  The  special  structure 
of  this  matrix  allows  us  to  compute  the  motion  and  structure  parameters  very  easily. 


2.  Preliminaries 

We  first  recall  some  details  about  perspective  projection,  the  motion  field,  the  brightness 
change  constraint  equation,  rigid  body  motion  and  planar  surfaces.  This  we  do  using 
vector  notation  in  order  to  keep  the  resulting  equations  as  compact  as  possible. 

2.1.  Perspective  Projection 

Let  the  center  of  projection  be  at  the  origin  of  a  Cartesian  coordinate  system.  Without 
loss  of  generality  we  assume  that  the  effective  focal  length  is  unity.  The  image  is  formed 
on  the  plane  2=1,  parallel  to  the  xy-plane,  that  is,  the  optical  axis  lies  along  the  2-axis. 
Let  R  be  a  point  in  the  scene.  Its  projection  in  the  image  is  r,  where 


The  2-component  of  r  is  clearly  equal  to  one,  that  is  r  •  z  =  1. 

2.2.  Motion  Field  and  Optical  Flow 

The  motion  field  is  the  vector  field  induced  in  the  image  plane  by  the  relative  motion  of 
the  observer  with  respect  to  the  environment.  The  optical  flow  is  the  apparent  motion  of 
brightness  patterns.  Under  favourable  circumstances  the  optical  flow  is  identical  to  the 
motion  field  (moving  shadows  or  uniform  objects  in  motion  could  create  descrepencies 
between  the  motion  field  and  the  optical  flow.  Here,  we  assume  that  the  motion  flow  field 
and  the  optical  flow  are  the  same).  The  velocity  of  the  image  r  of  a  point  R  is  given  by 

—  =  d  1  R. 
dt  dt  R  •  z 

For  convenience,  we  introduce  the  notation  rt  and  Rf  for  the  time  derivatives  of  r  and 
R,  respectively.  We  then  have 


rt 


1 

R  •  z 


R  t- 


l 

<r  • 


(R,  •  z)R, 


»  1  m  t,m\ 
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which  can  also  be  written  in  the  compact  form 


rt 


1 


(z  x  (R*  x  R)), 


since  a  x  (b  x  c)  =  (c  •  a)b  -  (a  •  b)c.  The  vector  r<  lies  in  the  image  plane,  and  so 
(r*  •  z)  =  0.  Further,  r*  =  0,  if  R*  ||  R,  as  expected. 

Finally,  noting  that  R  =  (R  •  z)r,  we  get 


r  t 


1 

R  ■  z 


(z  x  (R<  x  r)) . 


2.3.  Rigid  Body  Motion 

In  the  case  of  the  observer  moving  relative  to  a  rigid  environment  with  translational 
velocity  t  and  rotational  velocity  u,  we  find  that  the  motion  of  a  point  in  the  environment 
relative  to  the  observer  is  given  by 

R<  =  —  u  x  R  —  t. 

Since  R  =  (R  •  z)r,  we  can  write  this  as 

R*  =  -(R  •  z)  u>  x  r  —  t. 

Substituting  for  Rf  in  the  formula  derived  above  for  r t,  we  obtain 

n  =  -(zx(rx(rxw  -  — yrt))) . 

It  is  important  to  remember  that  there  is  an  inherent  ambiguity  here,  since  the  same 
motion  field  results  when  distance  and  the  translational  velocity  are  multiplied  by  an 
arbitrary  constant.  This  can  be  seen  easily  from  the  above  equation  since  the  same 
image  plane  velocity  is  obtained  if  one  multiplies  both  R  and  t  by  some  constant. 


2.4.  Brightness  Change  Equation 

The  brightness  of  the  image  of  a  particular  patch  of  a  surface  depends  on  many  factors. 
It  may  for  example  vary  with  the  orientation  of  the  patch.  In  many  cases,  however,  it 
remains  at  least  approximately  constant  as  the  surface  moves  in  the  environment.  If  we 
assume  that  the  image  brightness  of  a  patch  remains  constant,  we  have 


or 

dE  dr 

dr  '  dt  +  dt  ~  ’ 

where  dE/dr  =  (dE/dx,dE/dy,  0)T  is  the  image  brightness  gradient.  It  is  convenient 
to  use  the  notation  Er  for  this  quantity  and  E%  for  the  time  derivative  of  the  brightness. 
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Then  we  can  write  the  brightness  change  equation  in  the  simple  form 

Et  •  vt  +  Et  =  0. 

Substituting  for  r*  we  get 

Et  -  E,  •  (z  x  (r  x  (r  x  w  -  — 1  ,  t)))  =  0. 

Now 

£,-(ix(rxt))  =  ( Er  x  z)  •  (r  x  t)  =  (( Et  x  z)  x  r)  •  t, 
and  by  similar  reasoning 

E,  •  (z  x  (r  x  (r  x  w)))  =  (((EP  x  z)  x  r)  x  r)  •  w, 

so  we  have 

Et  -  ({[E,  x  z)  x  r)  x  r)  •  w  +  — T((Er  x  z)  x  r)  ■  t  =  0. 

1C  *  z 

To  make  this  constraint  equation  more  compact,  let  us  define  c  =  Et,  8  =  (Er  x  z)  xr, 
and  v  =  — s  x  r;  then,  finally, 


C  +  V-  W+  — — rS  •  t  =  0. 

R  •  z 

This  is  the  brightness  change  equation  in  the  case  of  rigid  body  motion. 

2.5.  Planar  Surface 

A  particularly  impoverished  scene  is  one  consisting  of  a  single  planar  surface.  The  equa¬ 
tion  for  such  a  surface  is 

R  •  n  =  1, 

where  n/  |n|  is  a  unit  normal  to  the  plane,  and  1/  |n|  is  the  perpendicular  distance  of  the 
plane  from  the  origin.  Since  R  =  (R  •  z)r,  we  can  write  this  as 

1 


so  the  constraint  equation  becomes 

c  -I-  v  •  u;  +  (r  •  n)(s  ■  t)  =0. 

This  is  the  brightness  change  equation  for  a  planar  surface.  Note  again  the  inherent 
ambiguity  in  the  constraint  equation.  It  is  satisfied  equally  well  by  two  planes  with  the 
same  orientation  but  at  different  distances  provided  that  the  translational  velocities  are 
in  the  same  proportions. 
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2.6.  Essential  Parameters  for  Planar  Surfaces 
The  brightness  change  equation  can  be  written  as: 

c  +  (r  x  s)  •  lj  +  (r  •  n)(s  •  t)  =  0. 

Now,  using  the  identity  (r  x  s)  •  u  =  (s  x  w)  •  r,  we  get: 

c  +  (s  x  u)  •  r  +  (r  •  n)(s  •  t)  =  0. 

Let  us  define: 

^  0  —  W3  +U>2  ^ 

W  =  1  +o>3  0  —u>i  I  , 

V  -u>2  +u>i  0  J 

then,  (s  x  w)  =  Ws,  in  which  case,  we  arrive  at: 

c  +  rTWs  +  rTntTs  =  0, 

and  after  collecting  terms: 

c  +  rT(W  +  ntr)s  =  0. 

Now  let: 

( Pi  P2  P3] 

P  =  I  P4  Ps  Pe  I  =  W  +  ntr, 

V  P7  P8  P9  ' 

We  will  refer  to  pt-  ,  i  =  1,2,  ...,9  as  the  essential  parameters  (in  agreement  with  Tsai  and 
Huang  [20])  since  these  parameters  contain  all  the  information  about  the  planar  surface 
orientation  and  motion  parameters.  Substituting  P  into  the  brightness  change  equation, 
we  get: 

c  +  rrPs  =  0. 

This  is  a  linear  constraint  equation  in  terms  of  the  elements  of  P,  and  can  be  used  to 
solve  for  these  parameters.  We  will  show  how  the  special  structure  of  P  can  be  exploited 
to  recover  the  motion  and  structure  parameters  very  easily. 

3.  Recovering  Essential  Parameters 

The  nine  essential  parameters  satisfy  the  following  constraint  equation: 

c  +  r^Ps  =  0. 

or  in  terms  of  p,  a  vector  of  length  9  whose  i-th  element  is  Pi : 

c  +  aTp  =  0, 

where: 

a  =  (riSi,riS2,riS3,r2Si,r2S2,’'2S3,r3Si,r3S2,r3S3):r. 

We  can  compute  them  using  image  brightness  E(x,y,t),  and  its  spatial  and  time  deriva¬ 
tives,  Er  and  Et ,  over  some  region  I  in  the  image  plane.  Since  there  are  only  eight  motion 
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and  surface  parameters  to  recover  (There  are  three  components  of  each  of  u> ,  t,  and  n. 
However,  as  mentioned  earlier,  the  translational  velocity  and  the  surface  normal  can  be 
recovered  up  to  a  scale  factor.),  only  eight  of  the  p,’s  are  independent.  This  implies  that 
we  can  arbitrarily  fix  one  of  the  essential  parameters,  and  compute  the  remaining  eight 
using  eight  independent  points  (At  each  point,  we  get  one  constraint  and  we  have  eight 
unknowns  to  recover). 

Let  p'  =  (pj ,  P2 ,  ...,Pg,0)  denote  the  solution  obtained  by  setting  the  last  element 
equal  to  zero.  If  we  define: 

p'  =  (p'i.P2>->P8)>  and  a  =  (riSi,r1s2,riS3,r2Si,r2S2,»,2«3,»'3Si,r3S2):r, 
then  the  above  constraint  equation  reduces  to: 

arp'  +  c  =  0. 

Using  eight  non-colinear  points  (this  is  necessary  to  guarantee  that  the  resulting  equations 
are  independent) ,  we  can  solve  the  following  linear  matrix  equation: 

Ap'  +  c  =  0, 


where: 


A  («*1 , 3-2 1  •••) as)  »  C  —  (ci,C2,...,Cg). 

The  solution  of  the  above  equation  is: 

p'  =  — A-1c. 


Image  brightness  values  are  distorted  with  sensor  noise  and  quantization  error.  These 
inaccuracies  are  further  accentuated  by  methods  used  for  estimating  the  brightness  gra¬ 
dient.  Thus  it  is  not  advisable  to  base  a  method  on  measurements  at  just  a  few  points. 
Instead  we  propose  to  minimize  the  error  in  the  brightness  constraint  equation  over  the 
whole  region  /  in  the  image  plane.  So  we  choose  the  essential  parameters  that  minimize: 


J J  (aTp'  +  c)2  dxdy. 


The  solution,  in  this  case,  is  given  by: 

p '  ~  ~(f f  aar  dxdy)~1(J J  cadxdy). 


Note  that,  in  general,  the  true  pg  is  nonzero.  We  can  show  that  the  solution  obtained 
through  the  assumption  that  P9  =  0,  p',  and  the  true  solution  (denoted  by  p)  are  related 
by  the  equation: 

P  =  p'  +  P9U,  u  =  (l,0,0,0,l,0,0,0,l)r. 

The  proof  goes  as  follows.  Since  8  =  [E,  x  z)  x  r,  then  rrs  =  r  •  ((£?P  x  z)  x  r)  =0.  For 
any  arbitrary  constant  /,  such  that  L  =  11  (I  is  the  identity  matrix),  we  have: 


rrLs  =  0. 


.•  -W-  /  V.V. 


.*  .v  •*  ^ 


.*  «*■ 
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If  we  add  this  to  our  costraint  equation: 

c  +  rT(W  +  ntr  +  L)s  =  0. 

It  is  immediately  apparent  that  any  P;  of  the  form: 

( p'i  P2  P'3\ 

P'  =  p'4  Ps  Pe  =  W  +  ntr  +  L 

V  P7  P8  Pg  / 

will  also  satisfy  our  constraint  equation.  Therefore,  the  two  solutions  for  the  essential 
parameters  are  related  by: 

P  =  P'  -  L, 

or  in  terms  of  p  and  p': 

p  =  p'-/u,  u  =  (1,0, 0,0, 1,0, 0,0,1), 

for  some  constant  /.It  follows  that  j>9  =  p9  —  l.  Now  if  p9  =  0,  then  l  =  — pg,  and  so: 

P  =  p'  +  Pgu. 

In  term  of  matrices  P,  and  P': 

P  =  P'  +  P9l. 

The  procedure  for  determining  pg,  and  subsequently  P  is  presented  in  the  appendix.  For 
now,  we  assume  that  P  is  known  (note  that  we  have,  so  far,  shown  how  to  compute  P' 
and  not  P). 


4.  Recovering  Motion  and  Structure 

We  now  show  how  to  compute  the  parameters  of  the  translational  motion  and  the  plane 
orientation  from  the  essential  parameters.  Once  these  are  known,  the  rotational  param¬ 
eters  are  determined  from: 

W  =  P  -  ntr 

Since  W  is  skew-symmetric,  and  from  the  definition  of  P: 

P*  =  P  +  Pr  =  {W  +  WT)  +  (tnr  +  ntr)  =  (tnr  +  ntr). 

Consider  the  special  case  where  t  ||  n  (n  =  kt,  for  some  constant  fc).  In  this  case, 

P‘  =  2kttT 

Computing  t  from  the  above  equation,  and  subsequently  n  is  a  simple  practice  of  algebra. 
Therefore,  from  now,  we  assume  that  t  x  n  /  0. 

Since  the  translational  vector  and  the  surface  normal  can  be  recovered  up  to  a  scale 
factor,  we  can,  without  loss  of  generality,  let  tTt  =  1.  We  will  present  our  results  in 
terms  of  the  normalized  n  and  t,  that  is: 

n  -  t 

n  =  r-r,  t  =  —  =  t. 

Inl  |t| 
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Note  that: 

TraceP*  —  (nTt  +  trn)  =  2(n  •  t)  =  2  jn|  (n  •  t), 
and  therefore,  we  can  always  copmute  n  from: 

■  I-  .I  TraceP* 

n  =  |n|n'  |n|  =  ^lT 

We  also  have  to  normalize  P*.  From  now,  P*  denotes  the  normalized  matrix,  so  that: 

P*  =  (tnT  -f  ntr). 

In  order  to  present  the  solutions  for  ii  and  t,  it  is  necessary  to  express  the  eigenvalue 
decomposition  of  P*  in  terms  of  these  vectors.  We  will  do  so  in  the  next  lemma. 

Lemma  1:  Let  P*  =  UAUr  be  the  eigenvalue  decomposition  of  P*  =  (tnr  +  ntr), 
and  let  r  =  ^-TraceP*.  Then: 

fr~l  0  0  \ 

A  =  0  0  0  , 

V  0  0  t  + 

JT_  t-n  txn  t  +  n 

L\/2(l-r))  jtxhj  V 2(1  +  t))J 

Proof:  (thT  +  ntT)t  x  n  =  t(n  •  (t  x  n))  +n(t  •  (t  x  n))  =  0  =  0  (t  x  n),  and  so  (t  x  n) 
is  the  eigenvector  corresponding  to  the  zero  eigenvalue  of  P*.  Since  it  is  symmetric,  P* 
has  3  orthogonal  eigenvectors.  The  other  two  eigenvectors  are,  therefore,  in  the  plane 
containing  t  and  n.  Let  u  =  at  -f  /?n  and  A  denote  an  eigenvector-eigenvalue  pair  for 
some  a  and  /3  (to  be  determined).  Then 

(tnr  +  ntr)(Qt  +  0n)  =  A(at  +  /?n), 

which  simplifies  to 

[a(t  •  n)  +  0(n  •  h)]t  +  (a(t  •  t)  +  P(t  ■  n)]ii  =  Aat  -(■  A/?n. 


Since  (t  •  n)  =  r,  we  have: 


The  solutions  for  A  are  given  by: 


A  =  r  ±  1. 


Substituting  for  A  into  the  earlier  equations,  we  get: 
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Using  these  into  the  equation  for  u,  and  normalizing  the  results  yield: 

t±n 

u  =  —  ,  . 

V 2(1  ±  t) 

Note  that  since  [r|  <  1,  Ai  <  A2  =  0  <  A3.  I 

We  can  now  detrmine  t  and  n.  Let  u,  denote  the  i-th  column  of  U.  From  lemma  1: 

t  -  h  t  +  n 


Ui  = 


V  2(1  —  t)  ’  U3  v/2(l  +  r)' 
From  the  expressions  for  the  eigenvectors,  it  follows  that: 

t  +  n  =  \/2(l  +  r)u3,  t  —  fi  =  \/  2(1  —  r)ui. 
Solving  these  for  t  and  n: 


,/  1  +  T 

ll  —  T  ~  ll+T  1 

n-V  2  u3  > 

/  2  Ul»  1  ~  V  2  U3  +  V 

Since  the  choice  of  the  signs  of  the  eigenvectors  are  arbitrary,  we  should  repeat  the  above 
procedure  with  the  sign  of  either  Ui  or  U3  reversed: 


Ur  =  -- 


t  -  n 


V  2(1  —  r)  ’  U3  \J  2(1  +  r) 

In  this  case,  the  second  solution  is  obtained  from  the  following  equations: 


t  +  n 


fl  +  T  ll~T  -  ll  +  T  /  1  —  T 

n=v— 1 «3  +  V— 1 “I.  t  =  \J  ~UZ-\J  — m, 

Note  that  this  is  the  dual  solution  for  the  planar  surfaces.  The  special  case  of  t  ||  n 
corresponds  to  when  the  matrix  P*  has  multiple  eigenvalues.  Then,  either: 

1)  r  =  1,  for  which  Aj  =  A2  =  0.  Then  the  two  solutions  merge  to  the  single  one: 

n  =  t  =  u3, 


or: 

2)  r  =  -1,  so  that  A2  =  A3  =  0.  Here,  the  unique  solution  is  given  by: 

n  =  -t  =  U|. 

As  mentioned  earlier,  we  can  determine  n  through  proper  scaling  of  n,  and  solve  for 
the  rotation  parameters  by  substituting  the  solutions  for  n  and  t  into  the  equation: 


W  =  P  -  ntr. 
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Even  though  we  gave  a  complete  and  compact  proof  of  the  dual  solution  in  an  earlier 
paper[l5],  it  is  intriguing  to  confirm  those  results  with  our  closed  form  solution.  We 
showed  that  the  two  solutions  are  related  by: 


n*  =  kt,  V  —  l/kn,  u*  =  u>  +  n  x  t, 

where  k  is  arbitrarily  chosen  to  scale  the  translation  vector  and  plane  normal.  The  two 
solutions  given  in  lemma  2  for  n  and  t  already  satisfy  the  duality  relationship  given 
above.  It  remains  to  show  the  same  for  the  rotation  parameters.  We  only  have  to  show 
that: 

W*  +  n*t*r  =  P, 


where 


or 


/  0  nit2-n2ti  nit3-niti\ 

W*  =  W  +  i  n2*i  -  ni*2  0  n2*3  -  723*2  .  , 

V  713*i  -  711*3  713*2  -  712*3  0  J 


W*  =  W  +  ntr  -  tnT. 


Substituting  for  n*,  t*,  and  W*  into  the  earlier  equation,  and  simplifying  the  results,  we 
get: 

W  +  ntr  =  P. 


5.  Summary 

In  this  paper,  we  presented  a  closed  form  solution  for  recovering  the  motion  of  an 
observer  with  respect  to  a  planar  surface  without  having  to  compute  the  optical  flow 
as  an  intermediate  step.  We  only  need  the  image  intensity  gradients  at  8  non-colinear 
points,  however,  in  general,  it  is  better  to  compute  these  gradients  in  a  larger  portion 
of  the  image  to  reduce  the  noise  effects.  We  first  employed  a  constraint  equation  we 
developed  for  planar  surfaces  to  compute  9  intermediate  parameters,  the  elements  of  a 
3x3  matrix.  We  referred  to  them  as  essential  parameters.  The  special  structure  of  this 
matrix  allows  us  to  compute  the  motion  and  plane  parameters  very  easily. 
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Appendix 

In  the  previous  sections,  we  showed  how  the  motion  paprameters  can  be  recovered 
once  the  essential  paprameters  are  known.  However,  the  brightness  change  constraint 
equation  allowed  us  to  determine  the  matrix  P'  (a  particular  solution  of  P  with  the  last 
element  set  to  zero).  We  showed  that  the  two  solution  are  related  by: 

P  =  P'  +  pgl. 

Here,  we  show  how  to  determine  pg,  and  consequently  P.  For  simplicity,  let  pg  =  1,  and 
let  P  =  TJAYt  denote  the  eigenvalue  decomposition  of  P,  where  (A  =  diag(Ai,  A2,  A3)) , 
then: 

P'  =  UAVr  -11  =  UAVt  -  /UVT. 

If  L  =  II  then: 

P'  =  UAVr  -  ULVr  =  U(A  -  L)Vr. 

Similarly,  if  P*  =  P  +  PT  =  UAU7,  denotes  the  eigenvalue  decomposition  of  P*,  then: 

P'*  =  P'  +  P'7’  =  P  +  PT  -  2L  =  P*  -  2L  =  U(A  -  2L)Ur. 

In  lemma  1,  we  showed  that  Ai  <  A2  =  0  <  A3  (and  when  t  jj  n,  we  get  two  zero 
eigenvalues) .  Therefore,  the  eigenvalues  of  P'*  can  be  arranged  in  the  form: 

Aj  —  2/  <  -21  <  A3  -  21. 

It  follows  that  l  =  -iA2(P'). 

So  in  summary,  we  assume  that  pg  =  0,  and  solve  for  the  essential  parameters  (ele¬ 
ments  of  P').  The  eigenvalue  decomposition  of  P'*  allows  us  to  determine  the  unknown 
shift  (pg  =  —  |A2),  and  then,  P  from: 

P  =  P'  -  ^A2(P')I 

Finally,  the  solution  of  the  motion  and  structure  parameters  are  determined  from  the 
equations  given  earlier  in  terms  of  the  trace  and  eigenvectors  of  P*  (note  that  the  eigen¬ 
vectors  of  P*  and  P'*  are  the  same). 
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